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Abstract 

In  this  paper,  we  consider  convergence  of  classical  high  order  Godnnov-type  schemes 
towards  entropy  solntions  for  scalar  conservation  laws.  It  is  well  known  that  sufficient  condi¬ 
tions  for  such  convergence  include  total  variation  boundedness  of  the  reconstruction  and  cell 
or  wavewise  entropy  inequalities.  We  prove  that  under  large  time  steps,  we  only  need  total 
variation  boundedness  of  the  reconstruction  to  guarantee  such  convergence.  We  discuss  high 
order  total  variation  bounded  reconstructions  to  fulhll  this  sufficient  condition  and  provide 
numerical  examples  on  one  dimensional  convex  conservation  laws  to  assess  the  performance  of 
such  large  time  step  Godnnov-type  methods.  To  demonstrate  the  generality  of  this  approach, 
we  also  prove  convergence  and  give  numerical  examples  for  a  large  time  step  Godunov-like 
scheme  involving  Sanders’  third  order  total  variation  diminishing  reconstruction  using  both 
cell  averages  and  point  values  at  cell  boundaries. 
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1  Introduction 


In  this  paper,  we  consider  the  Canchy  problem  of  the  scalar  conservation  law: 


Ut  +  V-fiu)=0,  in  M”^x[0,T] 
u{x,0)  =  uqIx),  in 


(1.1) 


where  uo{x)  G  BV,  the  space  of  fnnctions  with  bounded  variation.  We  do  not  consider 
boundary  conditions  and  assume  that  the  initial  condition  is  either  compactly  supported 
or  periodic.  In  most  of  our  presentation  we  will  concentrate  on  the  one  dimensional  case 
m  =  1,  and  will  point  out  similarities  and  differences  for  higher  dimensions  m  >  1.  The 
scheme  under  consideration  is  the  classical  Godunov-type  scheme,  consisting  of  the  following 
three  stages  to  evolve  the  cell  averages  from  time  level  n  (denoted  by  u'^)  to  time  level  n  -|-  1 
(denoted  by 

1.  Reconstruction:  obtain  a  high  order  piecewise  polynomial  reconstruction  u"'{x)  whose 
cell  averages  agree  with  the  given  cell  averages  u^.  We  denote  this  reconstruction 
operator  as  Re{vR). 

2.  Evolution:  evolve  u'^{x)  by  the  conservation  law  (1.1)  exactly  for  a  time  step  At,  to 
obtain  a  solution  u'^^^{x)  which  is  in  general  not  piecewise  polynomial  anymore.  We 
denote  this  evolution  operator  as 

3.  Averaging:  average  the  function  vA^^{x)  to  obtain  the  cell  averages  at  time  level 
n  -|-  1.  We  denote  this  averaging  operator  as 

The  Godunov-type  scheme,  using  the  notations  introduced  above,  can  be  described  ab¬ 
stractly  as 

=  A{uo)  .  . 

=  AoSAtoRe(u^),  n  =  0,l,--- ,iV-l.  ^  ^  ^ 

The  purpose  of  the  reconstruction  step  is  to  increase  the  formal  accuracy  of  the  scheme. 
It  is  usually  achieved  by  using  cell  averages  of  several  neighboring  cells  to  obtain  a  high  order 
reconstructed  polynomial.  In  order  to  be  more  specihc,  let  us  first  introduce  some  notations. 
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For  example,  in  the  one  dimensional  case,  the  mesh  is  given  by  /,•  =  ,  a;,- ,  il.  The  center 

of  the  cell  Ij  is  given  by  Xj  =  |  {^j-^  “*"^7+5)’  denoted  by 

Ax  =  x,,i  —  a;_i,  which  is  considered  to  be  a  constant  (uniform  mesh)  in  this  paper  for 
simplicity  of  notations.  The  cell  average  of  a  function  u{x)  in  the  cell  Ij  is  dehned  as 


%  =  A{u),  =  — 


'  X  .  1 


u{x)  dx 


A  second  order  reconstruction  would  be  a  piecewise  linear  function,  and  the  linear  recon¬ 
struction  pj  [x)  in  cell  Ij  could  be  the  unique  linear  polynomial  satisfying 


—  pj{x)dx  =  Uj_u  tT 


Ax 


Ax 


'  X  .  1 

•?--2 


Pj{x)  dx  =  Uj, 


namely,  Pj(x)  agrees  with  the  given  cell  averages  of  u{x)  over  the  two  cells  {Ij-i,  Ij},  which 
is  called  the  stencil  of  the  reconstruction  Pj{x).  Usually,  we  require  that  the  stencil  for  the 
reconstruction  in  cell  Ij  contains  at  least  Ij  itself.  Additional  cells  in  the  stencil  can  be 
chosen  according  to  symmetry  or  stability  (non-oscillation)  considerations.  For  example,  the 
stencil  for  the  linear  polynomial  Pj{x)  above  could  also  be  {Ij,  Ij+i},  which  would  of  course 
produce  a  different  reconstructed  linear  polynomial. 

The  reconstruction  is  called  (formally)  k-ih  order  accurate,  namely 


\Re{u){x)  —  u{x)\  <  CAx^,  x  G  Ij  (1.3) 

if  the  reconstructed  polynomial  in  Ij  is  (/c  — l)-th  degree  with  a  stencil  of  k  cells  over  which  the 
function  u{x)  is  smooth.  Here  and  below  C  (or  C)  will  denote  a  generic  constant  independent 
of  the  mesh  sizes  Ax  and  At,  which  may  have  different  values  in  different  locations.  For 
more  details  of  second  and  higher  order  reconstructions,  we  refer  to  [11,  18]. 

A  reconstruction  is  called  total  variation  bounded  (TVB)  if 


TV{Re{u))  <  (1  +  CAx)TV{u) 


(1.4) 


and  it  is  called  total  variation  diminishing  (TVD)  if  (1.4)  holds  with  U  =  0.  Here  the  total 
variation  of  a  piecewise  smooth  function  g{x)  is  dehned  as  usual 

TV{g)=  [  \g\x)\dx 
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and  the  total  variation  of  u  in  (1.4)  is  defined  with  u  considered  as  a  piecewise  constant 
function. 

Since  the  exact  evolution  operator  S' At  and  the  averaging  operator  A  do  not  increase  the 
total  variation  of  a  function,  the  total  variation  of  the  hnal  solution  of  the  Godunov-type 
scheme  at  a  hxed  time  T  =  NAt,  is  bounded  by  a  constant  C  independent  of  the  mesh  sizes 
Ax  and  At,  when  both  of  them  go  to  zero,  if  the  reconstruction  is  TVB.  This  constant  C  is 
the  total  variation  of  the  initial  data  if  the  reconstruction  is  TVD. 

By  compactness,  the  sequence  of  numerical  solutions  of  a  Godunov-type  scheme  with 
a  TVB  reconstruction,  when  the  mesh  sizes  Ax  and  At  both  go  to  zero,  has  a  conver¬ 
gent  subsequence  in  L^.  Since  the  scheme  is  conservative,  the  limit  of  such  a  convergent 
subsequence  is  a  weak  solution  to  (1.1).  However,  because  of  the  non-uniqueness  of  weak 
solutions,  we  cannot  conclude  the  convergence  of  the  original  sequence  of  numerical  solu¬ 
tions,  nor  can  we  claim  that  the  limit  of  a  convergent  subsequence  is  the  physically  relevant 
entropy  solution.  Typically,  we  would  need  an  additional  cell  entropy  inequality,  e.g.  [15], 
or  a  wavewise  entropy  inequality,  e.g.  [21],  to  imply  such  convergence.  Such  cell  or  wave- 
wise  entropy  inequalities  are  very  difficult  to  obtain  for  high  order  schemes,  especially  for 
non-convex  conservation  laws  for  which  such  entropy  inequalities  are  needed  for  all  entropies 
rather  than  just  for  one  of  them.  As  far  as  we  know,  within  self-similar  hnite  difference  or 
hnite  volume  framework,  such  cell  or  wavewise  entropy  inequalities  have  been  obtained  only 
for  second  order  MUSGL  type  schemes  for  scalar  one  dimensional  convex  conservation  laws 
[15,  21],  and  for  a  special  second  order  scheme  for  general  conservation  laws  which  evolves 
cell  averages  and  cell  slopes  independently  [1].  These  schemes  are  all  explicit  schemes  with 
the  usual  GFL  restriction  At  =  0(Ax).  It  seems  that  for  such  situation  a  cell  or  wavewise 
entropy  inequality  is  necessary  for  convergence  towards  entropy  solutions,  see  for  example 
[15]  in  which  it  is  proven  that  hnite  volume  schemes  (those  which  evolve  only  cell  averages 
and  rely  on  reconstruction  to  achieve  higher  order  accuracy)  which  satisfy  all  entropy  condi¬ 
tions  are  necessarily  only  hrst  order  accurate,  and  [20]  in  which  a  counter  example  is  shown 
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to  demonstrate  that  second  order  TVD  Godunov-type  scheme  can  converge  to  an  entropy 
violating  solution  for  non-convex  conservation  laws.  We  remark  that  cell  entropy  inequality 
for  the  square  entropy  can  be  obtained  in  more  general  situations  (arbitrary  order  of  accu¬ 
racy,  general  triangulation  in  multi-space  dimensions,  arbitrary  scalar  conservation  laws)  for 
the  semi-discrete  discontinuous  Galerkin  method  [7],  but  this  is  a  totally  different  class  of 
schemes  in  which  no  reconstruction  is  involved  and  the  complete  polynomial  within  a  cell  is 
evolved  in  time  to  obtain  high  order  accuracy. 

We  now  mention  a  few  papers  in  the  literature  which  directly  motivated  our  work  in 
this  paper.  Lions  and  Souganidis  [13]  proved  the  convergence  of  an  implicit  MUSGL  scheme 
to  the  entropy  solution  of  strictly  convex  scalar  one  dimensional  Hamilton- Jacobi  equation 
and  conservation  law  under  a  large  time  step  assumption,  namely  Ax  =  o{At).  For  one 
dimensional  convex  Hamilton-Jacobi  equations,  which  are  very  closely  related  to  conservation 
laws  (1.1),  Ferretti  [3]  proved  the  convergence  of  a  semi-Lagrangian  scheme  to  the  viscosity 
solution  (the  corresponding  solution  to  the  entropy  solution  for  conservation  laws),  again 
under  a  large  time  step  assumption  Ax  =  o{At).  This  convergence  holds  for  quite  general 
high  order  reconstructions,  e.g.  the  ENO  reconstruction  [6]  up  to  hfth  order  and  WENO 
reconstruction  [8]  up  to  ninth  order,  see  [3,  2].  The  message  in  [13,  3,  2]  seems  to  be  that 
it  is  easier  to  prove  convergence  for  large  time  step  {Ax  =  o{At))  schemes.  Of  course,  such 
choice  of  the  time  step  violates  the  traditional  GFL  condition  for  local  explicit  schemes.  One 
would  need  either  to  take  implicit  time  discretization,  as  in  [13],  or  semi-Lagrangian  time 
discretization,  namely  following  the  characteristics,  as  in  [3,  2]. 

In  this  paper  we  prove  that,  indeed,  the  large  time  step  assumption  Ax  =  o{At)  allows 
us  to  bypass  the  verihcation  of  cell  or  wavewise  entropy  conditions  and  prove  the  conver¬ 
gence  towards  the  entropy  solution  of  Godunov-type  schemes  solely  under  the  assumption  of 
a  TVB  reconstruction,  for  general  scalar  conservation  laws  (1.1).  Since  self-similar  TVB 
reconstruction  is  available  in  one  dimension  to  third  order  accuracy  [12],  we  then  have 
third  order  accurate  Godunov-type  schemes  in  one  dimension  which  can  be  proven  con- 
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vergent  towards  the  entropy  solution.  Of  course,  to  actually  implement  the  large  time  step 
Godunov-type  schemes  is  not  trivial.  We  discuss  an  effective  implementation,  similar  to  the 
semi-Lagrangian  scheme  in  [3,  2]  for  Hamilton- Jacobi  equations,  for  one  dimensional  scalar 
convex  conservations. 

The  paper  is  organized  as  follows.  Section  2  contains  our  main  convergence  result  and 
its  proof.  In  section  3,  the  TVB  property  of  the  numerical  solution  is  considered.  Examples 
of  high  order  TVB  reconstruction  are  given.  Section  4  discusses  an  effective  implementation 
of  the  large  time  step  Godunov-type  scheme  for  one  dimensional  scalar  convex  conservation 
laws  and  provides  numerical  results  verifying  the  performance  of  such  schemes.  In  order 
to  demonstrate  the  generality  of  this  approach,  we  consider  in  section  5  a  Godunov-like 
scheme  using  a  third  order  TVD  reconstruction  due  to  Sanders  [16],  which  uses  both  cell 
averages  and  point  values  at  cell  boundaries  for  the  reconstruction  and  does  not  exactly  fit 
our  large  time  step  Godunov-type  scheme  framework,  but  the  convergence  can  be  obtained 
along  similar  lines.  Goncluding  remarks  are  given  in  section  6. 

2  Convergence  of  large  time  step  Godunov-type  schemes 

We  consider  the  Godunov-type  scheme  (1.2)  for  solving  the  scalar  conservation  law  (1.1), 
under  the  following  three  assumptions: 

1.  The  reconstruction  operator  Re  satishes 

||h  -  Re{u)\\Li  <  CAxTV{u)  (2.1) 

for  any  piecewise  constant  function  -u; 

2.  The  total  variation  of  the  numerical  solution  is  uniformly  bounded 

TV{u^)  <  C  (2.2) 

for  all  mesh  sizes  Ax  and  At  satisfying  nAt  <  T ; 
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3.  Large  time  step  condition 


,  Ax 
lim  — —  =  0 
Ax^O  At 


(2.3) 


We  have  the  following  convergence  result  under  these  assumptions. 


Theorem  2.1.  The  numerical  solution  of  the  Godunov-type  scheme  (1.2)  for  solving  the 
scalar  conservation  law  (1.1)  converges  in  to  the  entropy  solution  of  (1.1)  when  the  mesh 
size  Ax  goes  to  zero,  if  the  three  assumptions  (2.1),  (2.2)  and  (2.3)  are  satished. 

Proof:  We  denote  the  exact  entropy  solution  of  (1.1)  at  time  level  by  u{-,t'^)  and  the 
difference  between  this  exact  solution  and  the  reconstructed  numerical  solution  =  Re(u^) 
at  time  level  by  e"". 

e^  =  \\ui;e)-Reiu^)\\Li.  (2.4) 

First,  we  notice  that 


iiA(«(-,n)-h-iUi  = 


< 


< 


\\A{u{;r))-AoS^toRe  (h’^-')IUi 
11^  {u{;r)-SMoRe{u^-^))  lUi 
||w(-,r)-^Aioi?e(h”-i)|Ui 
II^At  (M(-,t”-'))-^Atoi?e(h"-i)|Ui 

\\u{-,r-^)  -  Re{u^-^)\\L1  = 


(2.5) 


where  the  third  inequality  is  due  to  the  fact  that  the  averaging  operator  A  does  not  increase 
the  norm  of  a  function,  and  the  hfth  inequality  is  due  to  the  contraction  property  of 
the  exact  entropy  solution  operator  S^t  of  (1.1)  [10].  Also  for  a  bounded  variation  function 
tc,  it  is  easy  to  verify  that 

llw  -  A{w)\\li  <  CAxTV{w)  (2.6) 
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In  fact,  for  the  one  dimensional  case 


\w-A{w)\\li  =  [ 

j 

-  z/ 


\w[x 


—  Wj\dx  =  j 


w'iOd^ 


dx  <  ^ 


w{x)  —  w{xj)\dx 
w'{^)\d^dx 


Ax  E/  \w\0\d^  =  ^^TV{w) 


where  we  have  used  the  mean  value  theorem  to  conclude  that  the  cell  average  Wj  equals 
w{xj)  for  some  point  Xj  inside  Ij,  which  is  valid  if  w{x)  is  continuous  inside  Ij.  However, 
(2.6)  is  valid  also  when  w'{x)  contains  ^-functions,  as  long  as  it  is  of  bounded  variation.  The 
proof  for  the  multi- dimensional  case  is  similar. 

We  can  now  proceed  to  prove  the  desired  convergence 


e"  =  ||n(-,t”)-i?e(h’^)|Ui 

<  \\u{;n  -  M<,n)\\L^  +  ii^(«(-,n)  -hiui  +  \k  -  Re{u^)\\Li 

<  CAxTV{u{-,  r))  -f  +  CAxTV{u^) 

<  CAx  +  e^-^ 

<  ■■■ 

<  CnAx  +  e° 

""  /\.X 

=  CT  — - h  e°  — 0  as  Ax  0 

At 

where  for  the  third  inequality  we  have  used  (2.6),  (2.5)  and  the  assumption  (2.1),  for  the 
fourth  inequality  we  have  used  the  fact  that  the  exact  entropy  solution  is  of  bounded 

variation  and  the  assumption  (2.2),  and  for  the  last  limit  we  have  used  the  large  time  step 
assumption  (2.3).  ■ 


Next,  we  prove  that  the  assumption  (2.1)  is  satished  by  all  commonly  used  reconstruc¬ 
tions,  including  the  high  order  ENO  and  WENO  reconstructions  [6,  8].  As  before,  for 
simplicity  we  divide  the  computational  domains  into  uniform  cells  (in  one  dimension  these 


are  non-overlapping  intervals  of  the  same  length  Ax  covering  the  whole  computational  do¬ 
main).  We  require  the  reconstruction  to  be  a  piecewise  polynomial  function,  namely  it  is  a 
polynomial  over  each  cell  Ij.  The  reconstructed  polynomial  in  Ij  is  uniquely  determined  by 
requiring  that  it  has  the  same  cell  average  as  the  given  one  over  each  cell  in  a  stencil,  which 
is  a  set  of  hnitely  many  cells  (the  number  of  cells  in  this  set  is  fixed  and  does  not  depend 
on  the  mesh  size  Ax)  and  contains  the  given  cell  Ij.  For  such  reconstructions  we  have  the 
following  result. 


Lemma  2.2.  The  assumption  (2.1)  is  satished  by  any  piecewise  polynomial  reconstruction 
with  a  local  stencil  containing  hnitely  many  cells. 

Proof:  We  prove  this  lemma  for  the  one  dimensional  case.  The  proof  for  the  multi¬ 
dimensional  case  is  similar.  Assume  the  reconstruction  stencil  for  cell  Ij  is  Ij-i{j),  ■  ■  ■  ,  Ij+r{j)- 
Since  the  stencil  contains  Ij  and  has  only  hnitely  many  cells,  we  have  0  <  /(j),  r(j)  <  M 
for  a  hxed  integer  M.  The  reconstructed  polynomial  in  cell  Ij  can  be  written  as 


(2.7) 


where  Pj(0  are  polynomials  independent  of  Ax  and  the  location  Ij,  see  for  example  [18]  for 
the  explicit  expression  of  We  therefore  have  |pi(OI  ^  C  for  some  constant  C  when 

^  is  in  the  relevant  range  |.^|  <  Also,  by  consistency  (i.e.  the  reconstruction  is  exact  if 
Uj  =  1),  we  have 


P* 


1 


Hence 


p{x)  —  Uj 


E 


Pi 


i=-l{j) 


{Ujj^-i  Uj) 
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We  now  have 


|!i  -  =  2],  11“  “ 


E 


rU) 


Ax 


^j) 


dx 


< 


CM{2M  +  l)Aa;  E  \Uj+i  -  Uj 
j 

CAxTViu) 


Clearly,  Lemma  2.2  covers  both  the  linear  reconstructions,  for  which  the  reconstruction 
stencil  relative  to  Ij  is  hxed,  i.e.  /(j)  and  r(j)  do  not  depend  on  j,  as  well  as  the  ENO 
reconstruction  [6],  for  which  the  reconstruction  stencil  depends  on  the  local  smoothness  of 
the  data  u.  Since  the  WENO  reconstructions  [8]  are  linear  combinations  of  hxed  stencil 
reconstructions  with  bounded  weights  (linear  combination  coefficients,  which  depend  on  the 
local  smoothness  of  the  data  u),  Lemma  2.2  also  covers  WENO  reconstructions.  Notice  that 
Lemma  2.2  is  valid  for  reconstructions  of  arbitrarily  high  order  accuracy. 

With  Lemma  2.2  we  obtain  the  following  corollary. 

Corollary  2.3.  If  the  reconstruction  operator  Re  is  based  on  piecewise  polynomial  recon¬ 
struction  with  a  local  stencil  containing  hnitely  many  cells,  then  the  numerical  solution  of 
the  Godunov-type  scheme  (1.2)  for  solving  the  scalar  conservation  law  (1.1)  converges  in 
to  the  entropy  solution  of  (1.1)  when  the  mesh  size  Ax  goes  to  zero,  if  the  two  assumptions 
(2.2)  and  (2.3)  are  satished. 

Proof:  We  would  just  need  to  combine  the  results  in  Theorem  2.1  and  Lemma  2.2.  ■ 
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3  TVB  property  of  the  numerical  solution  and  exam¬ 
ples  of  the  TVB  reconstruction  operator 

The  result  of  the  previous  section  indicates  that,  apart  from  the  large  time  step  assumption 
(2.3),  which  causes  difficulties  in  the  actual  implementation  of  the  scheme  and  will  be  dis¬ 
cussed  in  more  detail  in  next  section,  we  would  only  need  to  have  a  total  variation  bounded 
numerical  solution  to  imply  convergence.  If  the  reconstruction  operator  Re  is  TVB,  i.e.  it 
satishes  (1.4),  then  we  easily  have 

TV(u^)  =  TV{AoSMoRe{u^-^)) 

<  TV{Re{u^-^)) 

<  {l+C/Ax)TV{u^-^) 

<  ■■■ 

<  {l+C/AxYTV{u^) 

<  {l+C/AtYTV{uo) 

<  e^^TViuo) 

for  any  n  and  At  such  that  nAt  <  T.  Here  we  have  used  the  fact  that  neither  the  averaging 
operator  A  nor  the  exact  entropy  solution  evolution  operator  S' At  increases  the  total  variation 
of  a  function  for  the  second  inequality,  and  the  assumption  (2.3)  for  the  sixth  inequality  (in 
fact,  here  we  do  not  need  the  large  time  step  assumption,  we  just  need  the  assumption 
that  Ax  <  CAt).  We  therefore  only  need  to  verify  the  total  variation  boundedness  of  the 
reconstruction  (1.4)  in  order  to  guarantee  convergence  towards  entropy  solutions  for  a  large 
time  step  Godunov-type  scheme. 

Next  we  discuss  two  examples  of  TVB  reconstructions  (1.4).  We  concentrate  our  attention 
on  self-similar  reconstructions,  namely  reconstructions  which  do  not  depend  explicitly  on  the 
mesh  size  Ax  (see  the  example  of  the  reconstruction  (2.7)  in  which  the  polynomials  Pi{^) 
do  not  depend  on  Ax  [18]).  In  fact,  almost  all  known  self-similar  TVB  reconstructions  are 
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actually  TVD.  It  is  very  difficult  to  obtain  self-similar  TVB  reconstructions  in  two  or  higher 
dimensions  [4],  hence  we  restrict  our  attention  in  one  dimension. 


1.  Second  order  MUSCL  reconstruction:  This  reconstruction  [19]  is  piecewise  linear  and  is 
second  order  accurate  except  at  extrema,  where  it  degenerates  to  hrst  order  accuracy. 
The  reconstructed  polynomial  in  cell  Ij  is  given  by 

p{x)  =  Uj  +  Sj  [x  —  Xj) 


where  the  slope  Sj  is  given  by 


Ujj^l  —  Uj  Uj  —  Uj-i 


Sj  =  m  .  .  ,  . 

^  '  Ax  Ax 


with  the  usual  dehnition  of  the  minmod  function  m  [5] 


sign(a)  +  sign(6)  . 

m{a,  b)  = - - - mm(|a|,  |o|). 


In  words,  the  slope  is  taken  as  the  smaller  one  (in  magnitude)  between  those  of  the  two 
reconstruction  linear  polynomials  based  on  the  stencils  {/j-i,  Ij}  and  {Jj,  T,+i},  nnless 
these  two  slopes  are  of  opposite  sign  (i.e.  we  are  at  an  extremnm  of  the  cell  average 
h),  in  which  case  the  slope  is  taken  as  zero.  This  reconstrnction  is  clearly  TVD. 


2.  Third  order  non- os  dilatory  central  reconstruction:  This  reconstruction  [12]  is  piece- 
wise  qnadratic  and  is  third  order  accnrate  except  at  extrema,  where  it  degenerates  to 
second  order  accnracy.  The  reconstructed  polynomial  in  cell  Ij  is  based  on  qj{x),  the 
reconstrncted  qnadratic  polynomial  over  the  stencil  Jj,  Jj+i},  with  modihcation 

to  satisfy  a  non-oscillatory  property  (no  new  extremnm)  and  a  maximnm  principle. 
The  hrst  step  is  to  obtain  the  non-oscillatory  property  and  is  given  by 


p]{x)  =  Uj  +  dj  {qj{x)  -  Uj) 


(3.1) 
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where  the  limiter  6i  is  dehned  as 


9,  =  ^ 


mm  I  — - 1  ,  , 

Mj  —  Uj  rrij  —  Uj 

M,_  1  —  Uj  m,,i  —  Ui 

mm  I  — - — , - 1  ,  , 

Mj  —  Ui  rrii  —  Ui 


where 


if  Uj-i  <  Uj  <  Uj+i 

if  Uj-i  >  Uj  >  Uj+i 

if  {uj+i  —  Uj)  (uj  —  Uj^i)  <  0 


Mj  =  maxgj(a;),  rrij  =  mingj(a;), 

x£lj 

Mj^i  =  max  (^{uj  +  %±i),  qj±i{Xj^i^  ,  =  min  (^{uj  +  Uj±i),  qj±i{Xj^i^ 

The  second  step  is  to  further  modify  p]  (x)  to  satisfy  a  maximum  principle.  If  we  denote 


M  =  maxMo(a;),  m  =  minMo(a;), 

xgij  x&R 


and 


M]  =  maxp](x), 

^  xeh 


m]  =  min  p](x), 


then  the  hnal  reconstructed  polynomial  is  given  by 


{p]{ 


X 


Pj{x)  =  { 


(1) 

{p]{x 
{p]{x 
{p]{x 


(2) 


-  Mj)  +  M, 

-  mj)  +  m, 

-  Mj)  +  M, 

-  m))  +  m, 


when  Mj  >  M  and  wj  >m 
when  Mj  <  M  and  m)  <  m 
when  Mj  >  M  and  m)  <  m 

when  Mj  >  M  and  mj  <  m 

when  Mj  <  M  and  m)  >m 

(3.2) 


with 


^(1)  _  M-Uj 


Mj  -  Uj '' 


^(2)  _  m  Uj 

m)  —  Uj 


It  is  proven  in  [12]  that  the  reconstruction  (3.2)  is  TVB.  Moreover,  we  have 


\\v]{.x)  -  MjIIli  =  OjWqj^x)  -  MjIIli, 


and 


\\pj{x)-Uj\\Li  =  e^9\\pi,(^x)-Uj\\Li,  i  =  l,2  or  \\pj{x)-Uj\\Li  =  \\p]{x)-Uj\\Li 
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with  0  <  <  1  and  0  <  <  1.  Therefore, 

WPjix)  -  hjIUi  <  \\p]{x)  -  UjWl^  <  \\qj{x)  -  UjWl^. 

Hence  clearly  assumption  (2.1)  is  also  satished  by  this  reconstruction,  using  Lemma 

2.2. 

Other  possible  TVB  reconstructions  in  one  dimension  include  the  very  high  order  TVD 
reconstructions  in  [14]  (however  these  reconstructions  degenerate  to  hrst  order  accuracy  at 
smooth  extrema),  or  the  high  order  TVB  reconstructions  in  [17]  (however  the  reconstructions 
are  no  longer  self-similar). 

4  Implementation  of  the  scheme  and  numerical  exam¬ 
ples  for  one  dimensional  scalar  convex  conservation 
laws 

In  this  section,  we  discuss  a  practical  implementation  of  the  large  time  step  Godunov-type 
scheme  discussed  in  the  previous  section,  for  one  dimensional  scalar  convex  conservation 
laws,  using  a  semi-Lagrangian  type  scheme  similar  to  that  used  in  [3,  2]  for  Hamilton-Jacobi 
equations. 

To  implement  one  time  step  of  the  Godunov-type  scheme,  assume  that  we  have  the  cell 

averages  u"'  available,  and  we  apply  a  reconstruction  operator  to  obtain  the  reconstruction 

u^{x)  =  Re{u"'){x),  which  is  a  piecewise  polynomial  function.  Our  job  would  be  to  evolve 

u^{x)  exactly  for  a  time  step  At,  to  obtain  a  solution  and  then  take  the  cell  averages 

of  this  solution  to  obtain  It  is  in  general  very  difficult  to  obtain  the  solution 

.^n+i(^),  especially  when  At  is  large,  however  we  do  not  need  the  details  of  this  solution 

.^n+i(^),  j^g  eell  averages  As  demonstrated  in  Figure  4.1,  if  we  draw  a  backward 

characteristic  line  T, ,  i  from  each  of  the  cell  boundaries  x, ,  i  at  time  level  back  to  time 

2  2 

level  t”  with  the  foot  located  at  t/j+i,  such  that  the  entropy  solution  of  (1.1)  is  a  constant 

u,- ,  1  along  this  backward  characteristic  line  T  ■ ,  i,  then  we  can  apply  the  divergence  theorem 
2  2  ‘ 
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over  the  trapezoid  region  hij  to  obtain  an  explicit  formula  for  as 

/y  •  1 

u'^{x)dx  - 
''i-h 

+  (a;/  - ij_.)  (4.1) 


Figure  4.1:  The  trapezoid  region  Vtj  in  the  x-t  plane. 


Now,  the  only  difficulty  left  for  the  implementation  of  the  Godunov-type  scheme  is  to 
hud  the  backward  characteristic  lines  Notice  that  there  might  be  multiple  forward 

characteristic  lines  from  the  reconstruction  in  the  previous  time  level  u^(x)  which  cross  at 
the  target  point  Even  though  only  one  of  those  characteristic  lines  (or  at  most 

two  if  the  solution  happens  to  be  discontinuous  at  should  be  correct  and  the 

remaining  forward  characteristic  lines  will  have  crossed  other  characteristic  lines  before  the 
time  and  should  be  ignored  as  they  are  no  longer  valid  after  such  crossing,  it  seems 
difficult  to  numerically  implement  the  Ending  of  this  correct  characteristic  line  with  low 
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computational  cost.  We  have  therefore  used  another  strategy  to  numerically  hnd  the  correct 

characteristic  line  F, ,  i.  First,  we  draw  forward  characteristic  lines  for  all  cell  L  from 

2 


and  (for  piecewise  polynomials  of  degree  at  least  2)  from  any  internal  extrema  of  u^{x)  inside 
the  interval  /*.  We  denote  the  left-most  and  right-most  intersections  of  these  characteristic 
lines  with  the  line  t  =  as  ^j_i  and  respectively,  and  dehne 

as  the  domain  of  influence  of  the  interval  for  the  time  level  .  Likewise,  we  draw  the 
two  forward  characteristic  lines  for  all  cell-boundary  points  from 

(v.)) .  (<.)) . 

and  denote  the  left  and  right  intersections  of  these  two  characteristic  lines  with  the  line 
t  =  as  respectively,  and  dehne 

as  the  domain  of  inhuence  of  the  cell-boundary  x^j^i  for  the  time  level  .  Clearly,  because 
of  the  convexity  of  f{u)  and  hence  the  monotonicity  of  f'{u),  the  sets 

C',+1  =  e  J,}  ,  i  =  I*  :  i  e  1  }  (4.4) 

contain  all  candidate  intervals  and  cell  boundary  points  at  time  level  which  will  contain 
forward  characteristic  lines  crossing  the  target  point  Furthermore,  we  can 

use  the  minimum  and  maximum  values  of  f'{u^)  to  estimate  the  largest  possible  region  of 

dependency  for  the  point  at  the  time  level  and  search  for  candidates  of 

the  sets  Cj^i  and  only  within  this  region  of  dependency  to  save  computational  cost. 
For  each  Jj  with  i  G  C',-,  i,  we  should  be  able  to  hnd  a  forward  characteristic  line  from 
for  some  point  x*  G  /*,  such  that  it  crosses  the  target  point  .  The 
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solution  u*  =  u'^{x*)  is  a  constant  along  this  line,  which  can  be  found  through  the  implicit 
characteristic  equation 

u*  =  ~  f  ■  (4-^) 

For  the  Burgers  equations  with  piecewise  polynomial  reconstruction  of  degree  no  larger  than 
4,  u*  can  be  obtained  analytically.  For  the  general  sitnation  a  Newton  iteration  might  be 
needed  to  hnd  u* .  Likewise,  for  each  with  i  G  we  can  hnd  a  forward  characteristic 

line  from  with  a  solution  value  u*  between 

crosses  the  target  point  The  solntion  u*  can  be  fonnd  by  the  the  characteristic 

eqnation  which  is  simpler  than  (4.5): 


/'(«*)  = 


^ 


i+4 


At 


(4.6) 


Among  all  the  {x*,u*}  and  {x^j^^,u*}  that  we  have  fonnd,  we  would  need  to  choose  the  foot 
of  the  correct  backward  characteristic  line  Uj^i  according  to  the  classical  Lax  formnla  [9] 


g[yj^)  =m.mg{y)  where  g{y)=  u'^{^)d^  +  AtL 


-  y 

At 


(4.7) 


where  the  lower  limit  —  cx)  in  the  integral  is  irrelevant  and  can  be  replaced  by  any  hxed  point, 
and 

L{v)  =  f  {fT\u)du 


is  the  Legendre  transform  of  f{u).  If  there  are  mnltiple  minimizers  to  g{y),  corresponding 
to  the  sitnation  that  the  solntion  is  discontinuous  at  j ,  then  we  will  just  take 

the  smallest  minimizer  as  yj+i,  corresponding  to  the  left  limit  valne  of  the  solntion  at  the 
discontinnity  (^Xj_^_i,t'^~^^'^.  It  is  here  that  we  use  the  property  that  the  conservation  law  has 
a  convex  flux  /(«),  as  the  Lax  formula  is  only  valid  for  convex  conservation  laws. 

We  now  snmmarize  onr  implementation  of  the  large  time  step  Godnnov-type  scheme  for 
a  scalar  one  dimensional  convex  conservation  law  (1.1),  from  to 


1.  Perform  a  reconstruction  to  obtain  a  piecewise  polynomial  u^{x). 
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2.  Find  the  domain  of  influence  Jj  for  each  cell  Jj  via  (4.2)  and  the  domain  of  influence 

for  each  cell  boundary  point  via  (4.3). 

3.  Find  the  set  C, ,  i  of  all  candidate  intervals  and  the  set  (F, ,  i  of  all  candidate  cell 

J'  2  2 

boundary  points  at  time  level  which  will  contain  forward  characteristic  lines  crossing 
the  target  point  j ,  via  (4.4). 

4.  For  each  Jj  with  i  G  C'j+i,  find  a  forward  characteristic  line  from  {x*,t'^)  with  the 

solution  value  u*,  via  (4.5).  Likewise,  for  each  with  i  G  Cj+i,  find  a  forward 
characteristic  line  from  with  a  solution  value  u*,  via  (4.6). 

5.  Use  the  Lax  formula  (4.7)  to  determine  the  correct  backward  characteristic  line  {i/j+i ,  Mj+i  } 
among  all  the  {x*,u*}  and  {x^_^_i,u*}  that  we  have  found  above. 

6.  Use  (4.1)  to  obtain  for  all  j. 

In  the  following  we  provide  numerical  results  of  the  large  time  step  Godunov-type  scheme 

2 

applied  to  the  Burgers  equation,  namely  (1.1)  with  f{u)  =  for  several  different  recon¬ 
structions  and  for  both  smooth  and  discontinuous  solutions.  We  take  At  =  in  all 

runs. 

Example  4.1.  Burgers  equation  with  the  initial  condition 

Uo{x)  =  sin(7ra;),  — 1  <  a;  <  1  (4.8) 

with  periodic  boundary  condition.  The  solution  at  t  =  0.2  is  still  smooth  and  we  record  the 
j^i  numerical  orders  of  accuracy  in  Table  4.1,  for  the  first  order  Godunov  method 

(no  reconstruction),  second  order  MUSGL  reconstruction,  and  third  order  non-oscillatory 
reconstruction  in  [12].  All  these  cases  correspond  to  TVB  reconstructions  and  hence  provable 
convergence  for  the  large  time  step  Godunov-type  scheme  according  to  the  discussions  in  the 
previous  two  sections.  We  can  see  clearly  from  Table  4.1  that  all  schemes  achieve  their 
designed  order  of  accuracy. 
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Table  4.1:  errors  and  numerical  orders  of  accuracy  for  different  TVD  reconstructions. 

Burgers  equation  with  initial  condition  (4.8).  t  =  0.2 


Number  of 

constant 

MUSGL 

Liu  and  Osher 

cells 

order 

order 

order 

34 

1.26E-2 

~ 

1.38E-3 

- 

1.41E-4 

- 

68 

6.30E-3 

1.00 

2.89E-4 

2.26 

2.75E-5 

2.35 

136 

3.80E-3 

0.72 

8.44E-5 

1.77 

4.38E-6 

2.65 

272 

2.00E-3 

0.92 

2.66E-5 

1.66 

5.85E-7 

2.90 

544 

1.06E-3 

0.91 

7.55E-6 

1.81 

7.65E-8 

2.93 

1088 

5.33E-4 

0.99 

1.76E-6 

2.09 

9.72E-9 

2.97 

Even  though  the  third  order  linear  central  reconstruction  (namely  the  stencil  for  recon¬ 
structing  the  polynomial  in  Ij  is  hxed  as  Jj,  Jj+i})  and  the  third  order  ENO  recon¬ 

struction  [6]  are  not  proven  TVB,  they  perform  well  in  this  smooth  test  case,  yielding  the 
expected  third  order  accuracy  for  the  large  time  step  Godunov-type  scheme,  see  Table  4.2. 


Table  4.2:  errors  and  numerical  orders  of  accuracy  for  third  order  linear  central  and  ENO 

reconstructions.  Burgers  equation  with  initial  condition  (4.8).  t  =  0.2 


Number  of 

central 

ENO 

cells 

Ll  g];];or 

order 

Ll  g];ror 

order 

34 

1.51E-4 

- 

1.37E-4 

- 

68 

2.65E-5 

2.51 

2.90E-5 

2.24 

136 

3.86E-6 

2.78 

5.05E-6 

2.52 

272 

5.21E-7 

2.88 

6.59E-7 

2.93 

544 

6.79E-8 

2.93 

9.63E-8 

2.77 

1088 

8.65E-9 

2.97 

1.21E-8 

2.98 

Example  4.2.  Burgers  equation  with  a  discontinuous  initial  condition 


Uo{x) 


1,  if  -  1  <  a;  <  0 

3,  if  0  <  a;  <  1 


(4.9) 


with  periodic  boundary  condition.  The  solutions  of  the  large  time  step  Godunov-type  scheme 
at  t  =  0.2  for  different  reconstructions  are  shown  in  Figure  4.2.  We  can  see  in  general  good 
agreement  between  the  numerical  and  the  exact  entropy  solutions.  The  first  order  Godunov 
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scheme  is  showing  some  stair-like  features  in  the  rarefaction  region,  which  is  well  known,  and 
the  schemes  based  on  linear  reconstructions,  namely  the  second  order  upwind  reconstruction 
(for  this  case  with  f'{u)  >  0  the  stencil  for  reconstructing  the  polynomial  in  Ij  is  hxed 
as  and  the  third  order  central  reconstruction  are  both  slightly  oscillatory.  The 

schemes  based  on  MUSCL,  Liu-Osher  non-oscillatory  and  ENO  reconstructions  are  all  non- 
oscillatory  as  expected. 

5  Convergence  of  the  large  time  step  Godunov-like 
scheme  based  on  Sanders’  TVD  reconstruction 

Sanders  [16]  proposed  a  third  order  accurate  TVD  (he  called  it  total  variation  non-expansive) 
scheme  by  evolving  both  the  cell  average  and  the  values  at  two  end  points  of  a  cell.  He  can 
then  obtain  a  third  order  TVD  scheme  instead  of  the  usual  Godunov-type  TVD  scheme  which 
degenerates  to  second  order  at  best.  This  scheme  does  not  fit  the  framework  of  Godunov- 
type  schemes  that  we  have  discussed  before,  since  it  evolves  not  just  the  cell  averages  but  also 
point  values  of  the  solution  at  cell  boundaries.  However,  the  solution  procedure  is  sufficiently 
similar  to  Godunov-type  schemes  and  therefore  we  will  call  it  a  Godunov-like  scheme.  The 
scheme  is  based  on  a  reconstruction  operator  Rs,  which  produces  a  piecewise  quadratic 
function.  After  a  suitable  initialization  to  obtain  a  piecewise  quadratic  initial  numerical 
solution  u^{x),  the  Godunov- type  scheme  based  on  Sanders’  reconstruction  operator  Rs  can 
be  summarized  as: 

1.  Evolution:  evolve  u'^{x)  by  the  conservation  law  (1.1)  exactly  for  a  time  step  At,  to 
obtain  a  solution  u^^^{x)  which  is  in  general  not  piecewise  polynomial  anymore.  Recall 
that  this  evolution  operator  is  denoted  as  S^t  iu^)- 

2.  Reconstruction:  obtain  a  piecewise  quadratic  polynomial  reconstruction  u^{x)  by  Sanders’ 
reconstruction.  We  denote  this  reconstruction  operator  as  =  Rs 
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Figure  4.2:  Burgers  equation  with  the  discontinuous  initial  condition  (4.9).  t  =  0.2.  80 
uniform  cells,  (a)  First  order  scheme  with  no  reconstruction;  (b)  second  order  MUSCL 
reconstruction;  (c)  second  order  upwind  reconstruction;  (d)  third  order  ENO  reconstruction; 
(e)  third  order  linear  central  reconstruction;  (f)  third  order  non-oscillatory  reconstruction  of 
Liu  and  Osher. 
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This  Godunov-like  scheme,  using  the  notations  introduced  above,  can  be  described  ab¬ 
stractly  as 

\u^  =  Rs{uq)  .  . 

=  RsoS^t{u^),  n  =  ,iV- 1.  ^  > 

We  will  hrst  discuss  the  TVD  reconstruction  operator  Rs,  then  prove  the  convergence  of 

the  large  time  step  Godnnov-like  scheme  based  on  this  reconstruction.  Finally  we  will  discuss 

the  implementation  of  this  scheme,  again  for  one  dimensional  scalar  convex  conservation  laws, 

and  provide  nnmerical  examples. 

5.1  Sanders’  TVD  reconstruction 


Sanders’  reconstrnction  Rs  maps  a  bonnded  variation  function  u{x)  to  a  piecewise  qnadratic 
polynomial,  satisfying  the  following  four  properties: 

1.  Agreement  of  the  cell  averages: 

/  Rs{u){x)dx  =  /  u{x)dx,  Vj  (5.2) 

Jij  Jij 

2.  TVD: 

TV{Rs(u))  <  TV{u)  (5.3) 

3.  Maximnm  (and  minimum)  principle: 

sup  Rs{u){x)  <  supu(a;),  inf  i?s(M)(a;)  >  infM(a;) 


4.  Third  order  accnracy: 


Rs{u){x)  =  u{x)  +  0(Aa;^) 


when  u  is  locally  smooth. 


Let  us  now  describe  the  details  of  the  reconstruction  operator  Rs.  First,  given  the  cell 
average  Uj  and  the  cell  boundary  valnes  u~^_i  and  u~  we  can  easily  construct  a  unique 
quadratic  polynomial  qj{x)  snch  that 


qj{x)  dx 


=  Uj, 
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In  fact,  if  we  denote 


^  2 


a  =  u:  I  —  tin,  h  =  u .  .  —  Un, 


(5.4) 


then  the  quadratic  polynomial  qj{x)  is  given  by 

X  ' 

qj{x)  =  q{a,  b){x)  +  uj  =  3(a  +  b)  (  ^  ^  )  +  (b  -  a) 


\  2 

X  —  X^\  ,,  ,  X  —  X 


Ax 


Ax 


^{a  +  b)+Uj  (5.5) 


We  denote  M  to  be  the  larger  one  (in  magnitude)  of  a  and  b  and  m  the  smaller  one,  P  =  ^, 

swp{u  —  Uj),  if  M  <  0 
if  M  >  0 


E  = 


inf(M  —  Uj), 

I  i 


and  E  =  Let 


T+  =  —E  ■ 


3(1 +  P) 

1  +  p  +  ’ 


=  ~2^  (P  +  3^)  -  3{E  -  p){3E  +  p)- 


and 


n  = 


Tr  = 


1, 

i/ 

1, 

i/ 

—  1  <  p  <  0  and 

a 

=  M 

min(r_,  1), 

i/ 

—  1  <  p  <  0  and 

b 

=  M 

min(r+,  1), 

i/ 

0  <  p  <  1 

1, 

i/ 

1, 

i/ 

—  1  <  p  <  0  and 

b 

=  M 

min(r_,  1), 

i/ 

—  1  <  p  <  0  and 

a 

=  M 

min(r+,  1), 

i/ 

0  <  p  <  1 

then  the  TVD  reconstruction  in  Ij  is  dehned  by 

Rs{u){x)  =  q{Tia,  Trb){x)  +  Uj 


(5.6) 


(5.7) 


with  q  defined  in  (5.5).  It  is  proved  in  [16]  that  the  reconstruction  Rs  is  third  order  accurate 
and  TVD. 


5.2  Convergence  of  the  Godunov-like  scheme  based  on  Sanders’ 
reconstruction 

In  this  subsection  we  prove  the  convergence  of  the  Godunov-like  scheme  (5.1)  based  on 
Sanders’  reconstruction. 
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Theorem  5.1.  The  numerical  solution  of  the  Godunov-like  scheme  (5.1)  with  Sanders’ 
reconstruction  (5.7)  for  solving  the  scalar  conservation  law  (1.1)  converges  in  to  the 
entropy  solution  of  (1.1)  when  the  mesh  size  Ax  goes  to  zero,  if  the  large  time  step  assumption 
(2.3)  is  satished. 

Proof:  We  again  denote  the  exact  entropy  solution  of  (1.1)  at  time  level  by  u{-,t'^)  and 
the  difference  between  this  exact  solution  and  the  numerical  solution  at  time  level 
by  e^: 

=  \\u{-,e)-u^L^.  (5.8) 

First,  the  total  variation  of  the  numerical  solution  is  bounded  (in  fact,  diminishing): 

TP(m^)  =  TV{Rso 

<  TV{u^-^) 

<  ■■■ 

<  TV{u^)  <  TViuo)  (5.9) 

where  we  have  used  the  fact  that  neither  the  reconstruction  operator  Rs  nor  the  exact 
entropy  solution  evolution  operator  S^t  increases  the  total  variation  of  a  function,  for  the 
second  inequality. 
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We  now  proceed  to  prove  the  desired  convergence 

e”  =  \\u{-,e) -u^Wli 

<  n  -  +  ||^Ai(w"-')  -  A  o  +  II A  o  SAtiu^-^)  -  M^IUi 

=  \\SAt{u{;e-^))  -  ^At(w”-')||Li  +  ||^Ai(w"-')  -  o  ^Ai(«”-')||Li 

+  ||/1  o  Rs  o  SAt{u^-^)  -  RS  o  ^At(M"-')|Ui 

<  \\u{-,r-^)  -  u^-^\\li  +  CAxTV{S^t{u^-^))  +  CAxTV{Rs  O  SAt{u^~^)) 

<  +  CAxTV{u^-^)  +  CAxTV{u^-^) 

<  e^-^  +  CAx 

<  ■■■ 

<  CnAx  +  e° 

""  /\  7T 

=  CT  — - h  e°  — 0  as  Ax  0 

At 

where  for  the  third  equality  we  have  used  the  cell  average  preserving  property  (5.2)  of  Rs] 
for  the  fourth  inequality  we  have  used  the  fact  that  the  exact  entropy  solution  operator  S' At 
is  an  contraction,  and  the  inequality  (2.6)  with  w  =  S'At('w”“^)  and  w  =  Rs  o  SAt{u"'~^) 
respectively;  for  the  hfth  inequality  we  have  used  the  fact  that  the  reconstruction  operator 
Rs  and  the  exact  entropy  solution  operator  S^t  are  both  TVD;  and  for  the  sixth  inequality 
we  have  used  the  fact  that  the  total  variation  of  the  numerical  solution  is  bounded  (5.9).  ■ 

5.3  Implementation  and  numerical  examples 

For  one  dimensional  scalar  convex  conservation  law  (1.1),  the  implementation  of  the  large 
time  step  Godunov-like  scheme  based  on  Sanders’  reconstruction  operator  Rs  is  similar  to 
that  for  the  large  time  step  Godunov-type  scheme  described  in  detail  in  section  4.  This  is 
because  the  procedure  described  in  section  4  provides  both  the  cell  average  and  the  cell 
boundary  values  Uj_i  and  Mj+i,  hence  a  and  b  in  (5.4)  are  available.  The  only  additional 
information  needed  is  the  relevant  extremum  value  E  in  (5.6),  which  is  in  general  difficult 
to  obtain.  We  follow  Sanders  [16]  and  take  instead  the  corresponding  extremum  value  E 
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at  the  previous  time  level  between  the  feet  of  the  characteristics  y,_i  and  y,-,  i.  This 
approximation  of  E  is  exact  if  there  is  no  discontinuity  within  the  trapezoid  Qj  (see  Figure 
4.1),  and  is  over-estimating  (larger  for  the  sup  and  smaller  for  the  inf)  if  there  is  at  least  one 
shock  within  Qj.  It  does  not  affect  the  accuracy  and  TVD  property  of  the  reconstruction 
Rs. 

In  the  numerical  experiments  we  take  At  =  Ax^'^  as  before.  We  test  the  large  time  step 
Godunov-like  scheme  based  on  both  a  unmodulated  third  order  reconstruction  {qj  as  given 
by  (5.5))  and  on  Sanders’  TVD  reconstruction  Rs.  We  again  hrst  test  the  smooth  solution 
as  in  Example  4.1.  The  errors  and  numerical  orders  of  accuracy  are  shown  in  Table  5.1. 
We  can  see  clearly  third  order  accuracy  for  the  scheme  based  on  both  reconstructions  from 
Table  5.1. 


Table  5.1:  errors  and  numerical  orders  of  accuracy  for  the  third  order  unmodulated  third 

order  reconstruction  (5.5)  and  Sanders’  TVD  reconstruction  Rs.  Burgers  equation  with 
initial  condition  (4.8).  t  =  0.2 


Number  of 

unmodulated  third  order 

Sanders 

cells 

error 

order 

error 

order 

40 

1.37E-5 

- 

1.37E-5 

- 

80 

2.26E-6 

2.59 

2.26E-6 

2.59 

160 

3.64E-7 

2.63 

3.64E-7 

2.63 

320 

5.12E-8 

2.83 

5.31E-8 

2.77 

640 

6.70E-9 

2.93 

6.78E-9 

2.97 

Next,  we  test  the  discontinuous  solution  as  in  Example  4.2.  Figure  5.1  displays  the 
numerical  solution.  We  can  see  again  good  agreement  between  the  numerical  and  the  exact 
entropy  solutions.  The  scheme  based  on  the  unmodulated  third  order  reconstruction  is 
slightly  oscillatory,  and  the  scheme  based  on  Sanders’  TVD  reconstruction  is  not  oscillatory, 
as  expected. 
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Figure  5.1:  Burgers  equation  with  the  discontinuous  initial  condition  (4.9).  t  =  0.2.  80 
uniform  cells.  Left:  unmodulated  third  order  reconstruction;  Right:  Sanders’  TVD  recon¬ 
struction. 

6  Concluding  remarks 

The  main  result  of  this  paper  is  the  proof  of  convergence  towards  entropy  solutions  of  large 
time  step  Godunov-type  schemes  with  total  variation  bounded  reconstructions  for  general 
scalar  conservation  laws.  It  is  somewhat  surprising  that  no  cell  or  wavewise  entropy  inequal¬ 
ities  are  needed  in  the  proof.  We  have  also  discussed  practical  implementation  of  such  large 
time  step  Godunov-type  schemes  for  scalar  one  dimensional  convex  conservation  laws  and 
have  provided  numerical  examples  to  demonstrate  the  performance  of  these  schemes.  Fu¬ 
ture  work  will  include  investigations  on  practical  implementation  for  more  general  equations, 
and  on  possible  elimination  of  the  total  variation  bounded  assumption  on  the  reconstruction 
(or  the  replacement  of  this  condition  by  a  weaker  condition  easier  to  satisfy  for  high  order 
reconstructions  in  multi-space  dimensions). 
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